library(dplyr)
library(tidyr)
library(ggplot2)
library(stringr)
library(lme4)
options(scipen=999)
df = read.csv('avg_tensor_by_roi_wide.csv',
colClasses = c('subject' = 'factor',
'site' = 'factor',
'visit' = 'factor')) %>%
rename(scan = visit)
outcomes = df %>%
select(where(is.numeric)) %>%
colnames
# replace ATROPOS measures with NA for select images (segmentation failed)
atropos_cols <- df %>%
select(contains('ATROPOS')) %>%
colnames()
df[(df$subject == '04001' & df$site == 'NIH' & df$scan == '01') |
(df$subject == '01003' & df$site == 'NIH' & df$scan == '01') |
(df$subject == '03002' & df$site == 'NIH' & df$scan == '02') |
(df$subject == '04003' & df$site == 'NIH' & df$scan == '01') |
(df$subject == '04001' & df$site == 'BWH' & df$scan == '02') |
(df$subject == '03001' & df$site == 'BWH' & df$scan == '01'), atropos_cols] <- NA
The purpose of the DTIMSCAMRAS study is to test for site effects on DTI metrics of different brain regions in the context of a traveling subjects design with a harmonized imaging protocol. Participants (n=11) with stable Multiple Sclerosis (MS) were scanned in 4 different locations: National Institute of Health (NIH), Brigham and Women’s Hospital (BWH), Johns Hopkins University (Hopkins), and The University of Pennsylvania (Penn). Two scans were collected per site visit, and participants were re-positioned between scans. The modalities that were collected include T1s (MPRAGE), FLAIR, and Diffusion-weighted images (DWI).
To preprocess imaging data, qsiprep
was used for bias correction, skull-stripping and co-registration of the
T1s and DWIs, and specialized denoising and artifact correction for the
DWIs. Then, DTIFIT was used to compute the following scalar
maps: * Fractional Anisotropy * Mean Diffusivity * Axial Diffusivity
* Radial Diffusivity
Concurrently, several segmentation pipelines were used on the T1s (+ FLAIRs in the case of MIMOSA) with the purpose of defining particular regions of interest (ROIs):
The segmentation results were used to average the scalar maps across all voxels belonging to each of the ROIs. The 36 resulting metrics make up the set of outcome variables which are the focus of the site effects test. Below, the outcome variables are listed.
data.frame(OUTCOME = outcomes) %>%
separate(OUTCOME, c("SEGMENTATION", "ROI", "SCALAR"), remove = FALSE)
Note the general
A permutation-based F-test was used to test for site effects in the absence of parametric assumptions. For each subject, the absolute difference between all possible pairs of intra-site and inter-site measures were computed. For instance, take the first row of the subject data.frame below, which has BWH as the “reference site”: the absolute differences between that row and the 6 non-BWH rows are computed for \(y\)——this process is repeated for all rows (and the results averaged) to arrive to the average inter-site difference for this subject; for the average intra-site differences the absolute difference of the 4 intrasite pairs are averaged.
df %>%
filter(subject == '01001') %>%
select(!where(is.numeric), ATROPOS_GM_AD) %>%
rename(y = ATROPOS_GM_AD)
These differences were averaged across all subjects resulting in the Mean Absolute Difference Between Sites (\(\overline{MAD_{B}}\)) and the Mean Absolute Difference Within Sites (\(\overline{MAD_{W}}\)). The test statistic is composed of their ratio:
\[ F = \frac{\overline{MAD_{B}}}{\overline{MAD_{W}}} \]
Key Findings: Permutation tests revealed significant site effects in that all 36 outcomes.
Subjects 03001 and 03002 didn’t complete imaging at all 4 sites. Subject 03001 completed imaging at Penn and BWH only and 03002 completed imaging at BWH, Hopkins, and the NIH only.
Subject 02001 is missing missing DTI data from their NIH visit and is therefore excluded from this analysis.
t_df <- df %>%
count(subject, site, .drop = FALSE) %>%
tidyr::pivot_wider(names_from = site, values_from = n) %>%
mutate(`Row Totals` = BWH + Hopkins + NIH + Penn)
t_df <- t_df %>%
summarize(across(where(is.integer), sum)) %>%
bind_rows(t_df, .) %>%
mutate(across(where(is.factor), as.character))
t_df[is.na(t_df)] <- "Column Totals"
t_df
Atropos segmentation failed for 6 images: 04001NIH01, 01003NIH01, 03002NIH02, 04003NIH01, 04001BWH02, and 03001BWH01.
atropos fail
#knitr::include_graphics(here::here(c('misc/atropos_fail.png', 'misc/atropos_success.png')))
#knitr::include_graphics("misc/atropos_fail.png")
#knitr::include_graphics("misc/atropos_success.png")
t_df <- df %>%
unite(sub, subject, site, scan, sep = '-') %>%
filter(!sub %in% c('04001-NIH-01', '01003-NIH-01', '03002-NIH-02', '04003-NIH-01', '04001-BWH-02', '03001-BWH-01')) %>%
separate(sub, c('subject', 'site', 'scan')) %>%
mutate_if(is.character, as.factor) %>%
count(subject, site, .drop = FALSE) %>%
tidyr::pivot_wider(names_from = site, values_from = n) %>%
mutate(`Row Totals` = BWH + Hopkins + NIH + Penn)
t_df <- t_df %>%
summarize(across(where(is.integer), sum)) %>%
bind_rows(t_df, .) %>%
mutate(across(where(is.factor), as.character))
t_df[is.na(t_df)] <- "Column Totals"
t_df
Densities are colored by site, while the black line represents the overall density aggregated across sites.
plot_densities <- function(var){
df %>%
ggplot(aes_string(x=var, color='site')) +
geom_density() +
stat_density(aes_string(x = var), geom = "line", color = "black") +
theme_bw() +
xlab(str_replace_all(var, "_", " "))
}
vars <-df %>%
select(ends_with('FA')) %>%
colnames()
dplots <- lapply(vars, plot_densities) %>%
setNames(vars)
for (name in names(dplots)){
cat("#####", str_replace_all(name, "_", " "), "\n")
print(dplots[[name]])
cat("\n\n")
}
plot_avg_tensor_values <- function(tensor){
df %>%
select(!is.numeric, ends_with(tensor)) %>%
pivot_longer(ends_with(tensor), names_to = 'seg', values_to = 'average') %>%
separate(seg, into = c('segmentation', 'roi', 'tensor')) %>%
unite(segmentation, segmentation, roi, sep = " ") %>%
ggplot(aes(x=site, y=average)) +
geom_boxplot(outlier.shape = NA) +
geom_jitter(aes(color = subject, shape=scan), alpha=0.6) +
#facet_grid(site~.) +
facet_grid(segmentation ~ .) +
coord_flip() +
ggtitle(sprintf("Mean %s in ROI by site", tensor)) +
theme_bw() +
# theme(axis.text.x = element_text(angle = 60, vjust = 1, hjust=1)) +
scale_x_discrete(expand = c(.00000001, .00000001)) +
xlab("") +
ylab("") +
theme(strip.text.y = element_text(angle = 0))
}
plot_avg_tensor_values('FA')
vars <-df %>%
select(ends_with('MD')) %>%
colnames()
dplots <- lapply(vars, plot_densities) %>%
setNames(vars)
for (name in names(dplots)){
cat("#####", str_replace_all(name, "_", " "), "\n")
print(dplots[[name]])
cat("\n\n")
}
plot_avg_tensor_values('MD')
vars <-df %>%
select(ends_with('RD')) %>%
colnames()
dplots <- lapply(vars, plot_densities) %>%
setNames(vars)
for (name in names(dplots)){
cat("#####", str_replace_all(name, "_", " "), "\n")
print(dplots[[name]])
cat("\n\n")
}
plot_avg_tensor_values('RD')
vars <-df %>%
select(ends_with('AD')) %>%
colnames()
dplots <- lapply(vars, plot_densities) %>%
setNames(vars)
for (name in names(dplots)){
cat("#####", str_replace_all(name, "_", " "), "\n")
print(dplots[[name]])
cat("\n\n")
}
plot_avg_tensor_values('AD')
pseudo-F Ratio statistics: Intersite Variability / Intrasite Variability
ratio_stats = readRDS("results/avg_tensor_by_roi_wide_ratio_stat.rds")
null_dist = readRDS("results/avg_tensor_by_roi_wide_null.rds")
ratio_stats
Distribution of permuted pseudo-F stats:
null_dist %>%
pivot_longer(everything()) %>%
ggplot(aes(x=name, y=value)) +
geom_violin(draw_quantiles = c(0.95)) +
geom_point(aes(x=name, y=value),
data = ratio_stats %>%
pivot_longer(everything())) +
coord_flip()
P-value of pseudo-F stat:
test_ratio_stat = function(ratio.stats, null.dists){
p.vals = c()
for (col in colnames(ratio.stats)){
ratio.stat = ratio.stats[, col]
null.dist = null.dists[, col]
percentile = ecdf(null.dist)
p.vals = c(p.vals, percentile(ratio.stat))
}
names(p.vals) = colnames(ratio.stats)
return(p.vals)
}
p = 1 - test_ratio_stat(ratio_stats, null_dist)
p_df <- data.frame(p_val = (10000*p + 1)/(1+10000))
p_df %>%
mutate(p_val_fdr = p.adjust(p_val, method = 'fdr'))
models <- outcomes %>%
lapply(function(x) lmer(reformulate("(1|subject) + (1|site:subject)", x), data=df)) %>%
setNames(outcomes)
icc_df <- models %>%
lapply(performance::icc, by_group = TRUE) %>%
bind_rows(.id = "outcome") %>%
as.data.frame()
icc_df %>%
filter(Group == 'subject') %>%
ggplot(aes(x=outcome, y=ICC)) +
geom_point() +
coord_flip() +
ggtitle("ICC of Subject Random Effect") +
ylim(0,1)
icc_df %>%
filter(Group == 'site:subject') %>%
ggplot(aes(x=outcome, y=ICC)) +
geom_point() +
coord_flip() +
ggtitle("ICC of Site:Subject Random Effect") +
ylim(0,1)
perf_mixed_model = function(model){
list(summary = summary(model),
perf = performance::icc(model, by_group = TRUE),
shapiro.resid = shapiro.test(residuals(model)),
shapiro.subject = shapiro.test(coef(model)$`subject`[,1]),
`shapiro.site:subject` = shapiro.test(coef(model)$`site:subject`[,1]),
qq.resid = car::qqPlot(residuals(model),dist="norm", main="Residuals"),
qq.ranef.subject = car::qqPlot(ranef(model)$subject[,1],dist="norm", main="Rand. Eff. subject"),
`qq.ranef.site:subject` = car::qqPlot(ranef(model)$`site:subject`[,1],dist="norm", main="Rand. Eff. site:subject"),
res.v.fit = ggplot(mapping = aes(y = resid(model), x = fitted(model))) + geom_abline(intercept = 0, slope = 0) + geom_point() + geom_smooth() + theme_bw()
)
}
perf_mixed_model(models$ATROPOS_WM_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_WM_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -507.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2705 -0.4084 -0.0074 0.4260 1.8737
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00004437 0.006661
## subject (Intercept) 0.00042895 0.020711
## Residual 0.00001039 0.003223
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.387130 0.006349 60.97
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.092
## subject | 0.887
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.97572, p-value = 0.1649
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.94623, p-value = 0.5962
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.9703, p-value = 0.3677
##
##
## $qq.resid
## 53 37
## 51 36
##
## $qq.ranef.subject
## [1] 11 1
##
## $`qq.ranef.site:subject`
## [1] 8 18
##
## $res.v.fit
perf_mixed_model(models$ATROPOS_GM_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_GM_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -602.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1815 -0.4246 0.0087 0.2757 4.3686
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000005326 0.002308
## subject (Intercept) 0.000038020 0.006166
## Residual 0.000005805 0.002409
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.160945 0.001919 83.87
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.108
## subject | 0.774
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.82302, p-value = 0.0000000553
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.98866, p-value = 0.9958
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.96679, p-value = 0.2835
##
##
## $qq.resid
## 51 52
## 49 50
##
## $qq.ranef.subject
## [1] 9 2
##
## $`qq.ranef.site:subject`
## [1] 8 18
##
## $res.v.fit
perf_mixed_model(models$FAST_WM_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_WM_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -557.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.23611 -0.48150 0.01147 0.42822 2.17374
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00002158 0.004646
## subject (Intercept) 0.00057907 0.024064
## Residual 0.00001383 0.003719
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.403259 0.007307 55.19
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.035
## subject | 0.942
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.97818, p-value = 0.1882
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93861, p-value = 0.5043
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.97583, p-value = 0.5383
##
##
## $qq.resid
## [1] 32 31
##
## $qq.ranef.subject
## [1] 1 11
##
## $`qq.ranef.site:subject`
## [1] 24 8
##
## $res.v.fit
perf_mixed_model(models$FAST_GM_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_GM_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -639.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4185 -0.3794 0.0231 0.3556 4.4690
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000006098 0.002469
## subject (Intercept) 0.000059127 0.007689
## Residual 0.000006814 0.002610
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.173461 0.002371 73.15
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.085
## subject | 0.821
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.84112, p-value = 0.00000008246
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.94055, p-value = 0.5269
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.98777, p-value = 0.9369
##
##
## $qq.resid
## [1] 51 52
##
## $qq.ranef.subject
## [1] 9 2
##
## $`qq.ranef.site:subject`
## [1] 8 18
##
## $res.v.fit
perf_mixed_model(models$JLF_WM_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_WM_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -570.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.15403 -0.62253 0.05912 0.47991 1.65535
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000050221 0.007087
## subject (Intercept) 0.000468884 0.021654
## Residual 0.000006638 0.002576
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.426441 0.006635 64.27
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.096
## subject | 0.892
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.99012, p-value = 0.8027
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.88581, p-value = 0.1232
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.95946, p-value = 0.1605
##
##
## $qq.resid
## [1] 53 54
##
## $qq.ranef.subject
## [1] 11 2
##
## $`qq.ranef.site:subject`
## [1] 8 24
##
## $res.v.fit
perf_mixed_model(models$JLF_GM_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_GM_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -671.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.32410 -0.42748 -0.05894 0.43779 2.37595
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000015087 0.003884
## subject (Intercept) 0.000041948 0.006477
## Residual 0.000002227 0.001492
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.149719 0.002058 72.75
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.255
## subject | 0.708
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.96707, p-value = 0.03697
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93611, p-value = 0.4759
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.96696, p-value = 0.2871
##
##
## $qq.resid
## [1] 78 77
##
## $qq.ranef.subject
## [1] 9 2
##
## $`qq.ranef.site:subject`
## [1] 28 27
##
## $res.v.fit
perf_mixed_model(models$FIRST_THALAMUS_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FIRST_THALAMUS_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -510.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.49646 -0.51505 0.04177 0.48047 1.81130
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00007717 0.008785
## subject (Intercept) 0.00036529 0.019113
## Residual 0.00002168 0.004656
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.309516 0.005958 51.95
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.166
## subject | 0.787
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98419, p-value = 0.4289
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.9135, p-value = 0.2682
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.96969, p-value = 0.3518
##
##
## $qq.resid
## [1] 52 70
##
## $qq.ranef.subject
## [1] 2 7
##
## $`qq.ranef.site:subject`
## [1] 24 26
##
## $res.v.fit
perf_mixed_model(models$JLF_THALAMUS_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_THALAMUS_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -464.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6797 -0.4081 0.0280 0.3874 3.4486
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00019846 0.014088
## subject (Intercept) 0.00023697 0.015394
## Residual 0.00003798 0.006163
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.319975 0.005212 61.4
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.419
## subject | 0.501
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.90137, p-value = 0.00001392
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.90554, p-value = 0.2158
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.96625, p-value = 0.2722
##
##
## $qq.resid
## [1] 51 52
##
## $qq.ranef.subject
## [1] 1 9
##
## $`qq.ranef.site:subject`
## [1] 10 8
##
## $res.v.fit
perf_mixed_model(models$MIMOSA_LESION_FA)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: MIMOSA_LESION_FA ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -385.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.7344 -0.3797 0.0021 0.3696 3.2023
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.0003639 0.01908
## subject (Intercept) 0.0013332 0.03651
## Residual 0.0001144 0.01070
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.3258 0.0115 28.34
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.201
## subject | 0.736
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.93547, p-value = 0.0005552
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.97161, p-value = 0.9023
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.87229, p-value = 0.0003283
##
##
## $qq.resid
## [1] 24 23
##
## $qq.ranef.subject
## [1] 4 5
##
## $`qq.ranef.site:subject`
## [1] 14 33
##
## $res.v.fit
perf_mixed_model(models$ATROPOS_WM_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_WM_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1549.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5582 -0.3564 0.0353 0.4075 2.8736
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000000000046271 0.000006802
## subject (Intercept) 0.000000000236783 0.000015388
## Residual 0.000000000004724 0.000002174
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000616523 0.000004775 129.1
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.161
## subject | 0.823
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.81595, p-value = 0.00000003449
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93972, p-value = 0.5172
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.95053, p-value = 0.0791
##
##
## $qq.resid
## 53 54
## 51 52
##
## $qq.ranef.subject
## [1] 1 3
##
## $`qq.ranef.site:subject`
## [1] 18 21
##
## $res.v.fit
perf_mixed_model(models$ATROPOS_GM_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_GM_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1503.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.65330 -0.40529 0.02173 0.35056 2.27089
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000007871 0.000008872
## subject (Intercept) 0.00000000020654 0.000014371
## Residual 0.00000000001163 0.000003410
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000754518 0.000004583 164.6
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.265
## subject | 0.696
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.95972, p-value = 0.0188
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.90487, p-value = 0.2118
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.97157, p-value = 0.403
##
##
## $qq.resid
## [1] 14 13
##
## $qq.ranef.subject
## [1] 8 1
##
## $`qq.ranef.site:subject`
## [1] 32 12
##
## $res.v.fit
perf_mixed_model(models$FAST_WM_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_WM_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1671.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5150 -0.4235 -0.0098 0.4002 2.8842
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000000000063511 0.000007969
## subject (Intercept) 0.000000000262843 0.000016212
## Residual 0.000000000005098 0.000002258
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000611772 0.000005061 120.9
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.192
## subject | 0.793
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.86709, p-value = 0.000000636
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93033, p-value = 0.4142
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.91937, p-value = 0.007356
##
##
## $qq.resid
## [1] 53 54
##
## $qq.ranef.subject
## [1] 1 3
##
## $`qq.ranef.site:subject`
## [1] 18 21
##
## $res.v.fit
perf_mixed_model(models$FAST_GM_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_GM_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1644.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.57814 -0.48669 -0.01082 0.39656 2.24634
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000006816 0.000008256
## subject (Intercept) 0.00000000016182 0.000012721
## Residual 0.00000000001044 0.000003231
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000736558 0.000004076 180.7
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.284
## subject | 0.673
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.96328, p-value = 0.02141
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.91875, p-value = 0.3084
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.98161, p-value = 0.7491
##
##
## $qq.resid
## [1] 14 13
##
## $qq.ranef.subject
## [1] 1 8
##
## $`qq.ranef.site:subject`
## [1] 18 22
##
## $res.v.fit
perf_mixed_model(models$JLF_WM_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_WM_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1655
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.10744 -0.47068 0.05638 0.40899 2.08380
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000002550 0.000005050
## subject (Intercept) 0.00000000028330 0.000016832
## Residual 0.00000000001291 0.000003593
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000649895 0.000005156 126
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.079
## subject | 0.881
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98922, p-value = 0.7456
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.91596, p-value = 0.2864
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.97438, p-value = 0.4894
##
##
## $qq.resid
## [1] 53 35
##
## $qq.ranef.subject
## [1] 1 3
##
## $`qq.ranef.site:subject`
## [1] 18 22
##
## $res.v.fit
perf_mixed_model(models$JLF_GM_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_GM_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1624.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.6447 -0.3965 0.0703 0.4914 1.8407
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000000000121841 0.000011038
## subject (Intercept) 0.000000000461916 0.000021492
## Residual 0.000000000008956 0.000002993
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000802476 0.000006729 119.3
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.206
## subject | 0.779
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.99124, p-value = 0.8673
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.92387, p-value = 0.3522
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.97873, p-value = 0.6423
##
##
## $qq.resid
## [1] 58 57
##
## $qq.ranef.subject
## [1] 8 1
##
## $`qq.ranef.site:subject`
## [1] 18 6
##
## $res.v.fit
perf_mixed_model(models$FIRST_THALAMUS_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FIRST_THALAMUS_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1585.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.21713 -0.55661 0.09266 0.56833 2.12975
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000006565 0.000008103
## subject (Intercept) 0.00000000021540 0.000014676
## Residual 0.00000000003787 0.000006154
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000664029 0.000004667 142.3
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.206
## subject | 0.675
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98821, p-value = 0.6786
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.92793, p-value = 0.3903
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.96794, p-value = 0.309
##
##
## $qq.resid
## [1] 21 35
##
## $qq.ranef.subject
## [1] 3 7
##
## $`qq.ranef.site:subject`
## [1] 22 24
##
## $res.v.fit
perf_mixed_model(models$JLF_THALAMUS_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_THALAMUS_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1574.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.07380 -0.40130 0.02501 0.52883 1.57687
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000009030 0.000009503
## subject (Intercept) 0.00000000030257 0.000017394
## Residual 0.00000000003792 0.000006158
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000651264 0.000005509 118.2
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.210
## subject | 0.702
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.97973, p-value = 0.2348
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93132, p-value = 0.4243
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.98705, p-value = 0.9206
##
##
## $qq.resid
## [1] 51 21
##
## $qq.ranef.subject
## [1] 1 3
##
## $`qq.ranef.site:subject`
## [1] 22 29
##
## $res.v.fit
perf_mixed_model(models$MIMOSA_LESION_MD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: MIMOSA_LESION_MD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1352.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.90417 -0.44509 -0.02298 0.38772 1.58995
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000000002872 0.00005359
## subject (Intercept) 0.000000003487 0.00005905
## Residual 0.000000000462 0.00002149
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.00091264 0.00001993 45.79
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.421
## subject | 0.511
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.99023, p-value = 0.8093
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.94465, p-value = 0.5766
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.9748, p-value = 0.5034
##
##
## $qq.resid
## [1] 24 19
##
## $qq.ranef.subject
## [1] 6 7
##
## $`qq.ranef.site:subject`
## [1] 33 14
##
## $res.v.fit
perf_mixed_model(models$ATROPOS_WM_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_WM_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1462.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6097 -0.2817 0.0546 0.2579 2.7148
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000018222 0.000013499
## subject (Intercept) 0.00000000007700 0.000008775
## Residual 0.00000000002234 0.000004727
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000891481 0.000003463 257.4
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.647
## subject | 0.273
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.81743, p-value = 0.00000003804
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93128, p-value = 0.4239
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.87728, p-value = 0.0004449
##
##
## $qq.resid
## 53 54
## 51 52
##
## $qq.ranef.subject
## [1] 1 3
##
## $`qq.ranef.site:subject`
## [1] 18 17
##
## $res.v.fit
perf_mixed_model(models$ATROPOS_GM_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_GM_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1488
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.39141 -0.38565 0.05608 0.39620 2.17972
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000013134 0.000011460
## subject (Intercept) 0.00000000016134 0.000012702
## Residual 0.00000000001253 0.000003539
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000876666 0.000004272 205.2
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.430
## subject | 0.529
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.97764, p-value = 0.2132
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93937, p-value = 0.5131
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.98065, p-value = 0.7137
##
##
## $qq.resid
## [1] 14 13
##
## $qq.ranef.subject
## [1] 1 8
##
## $`qq.ranef.site:subject`
## [1] 18 22
##
## $res.v.fit
perf_mixed_model(models$FAST_WM_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_WM_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1589.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5799 -0.3672 0.0311 0.3119 2.6620
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000016955 0.000013021
## subject (Intercept) 0.00000000007885 0.000008880
## Residual 0.00000000002294 0.000004789
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000898170 0.000003435 261.4
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.625
## subject | 0.291
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.87578, p-value = 0.000001328
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93498, p-value = 0.4634
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.88495, p-value = 0.0007178
##
##
## $qq.resid
## [1] 53 54
##
## $qq.ranef.subject
## [1] 3 1
##
## $`qq.ranef.site:subject`
## [1] 18 21
##
## $res.v.fit
perf_mixed_model(models$FAST_GM_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_GM_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1624.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1948 -0.4416 -0.0089 0.4843 1.9868
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000011830 0.00001088
## subject (Intercept) 0.00000000010126 0.00001006
## Residual 0.00000000001239 0.00000352
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000865994 0.000003523 245.8
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.510
## subject | 0.437
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.99014, p-value = 0.8036
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.90643, p-value = 0.2211
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.97625, p-value = 0.5527
##
##
## $qq.resid
## [1] 14 13
##
## $qq.ranef.subject
## [1] 1 8
##
## $`qq.ranef.site:subject`
## [1] 18 22
##
## $res.v.fit
perf_mixed_model(models$JLF_WM_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_WM_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1552.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.44605 -0.47717 0.00559 0.50950 1.50231
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000016048 0.000012668
## subject (Intercept) 0.00000000028930 0.000017009
## Residual 0.00000000004579 0.000006767
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000973298 0.000005573 174.7
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.324
## subject | 0.584
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98631, p-value = 0.5545
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.95076, p-value = 0.6537
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.93257, p-value = 0.01955
##
##
## $qq.resid
## [1] 53 21
##
## $qq.ranef.subject
## [1] 1 2
##
## $`qq.ranef.site:subject`
## [1] 18 24
##
## $res.v.fit
perf_mixed_model(models$JLF_GM_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_GM_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1590
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.72405 -0.41786 0.02606 0.49010 1.64221
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000020645 0.000014368
## subject (Intercept) 0.00000000037544 0.000019376
## Residual 0.00000000001495 0.000003866
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000918368 0.000006299 145.8
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.346
## subject | 0.629
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98905, p-value = 0.734
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.94762, p-value = 0.6137
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.95376, p-value = 0.1022
##
##
## $qq.resid
## [1] 53 63
##
## $qq.ranef.subject
## [1] 1 8
##
## $`qq.ranef.site:subject`
## [1] 18 20
##
## $res.v.fit
perf_mixed_model(models$FIRST_THALAMUS_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FIRST_THALAMUS_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1515.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.41831 -0.45695 0.06954 0.55154 1.63084
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000016471 0.000012834
## subject (Intercept) 0.00000000060854 0.000024669
## Residual 0.00000000008871 0.000009419
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000881751 0.000007795 113.1
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.191
## subject | 0.706
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.97775, p-value = 0.1769
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.95359, p-value = 0.69
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.96804, p-value = 0.3114
##
##
## $qq.resid
## [1] 21 53
##
## $qq.ranef.subject
## [1] 3 1
##
## $`qq.ranef.site:subject`
## [1] 24 14
##
## $res.v.fit
perf_mixed_model(models$JLF_THALAMUS_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_THALAMUS_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1520
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.52298 -0.38570 0.04951 0.55046 1.54271
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000015241 0.000012346
## subject (Intercept) 0.00000000064654 0.000025427
## Residual 0.00000000008218 0.000009065
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000873952 0.000007988 109.4
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.173
## subject | 0.734
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.96515, p-value = 0.02801
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.94975, p-value = 0.6407
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.96767, p-value = 0.303
##
##
## $qq.resid
## [1] 21 53
##
## $qq.ranef.subject
## [1] 3 2
##
## $`qq.ranef.site:subject`
## [1] 24 39
##
## $res.v.fit
perf_mixed_model(models$MIMOSA_LESION_AD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: MIMOSA_LESION_AD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1327.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.89760 -0.40883 -0.00999 0.31963 1.92297
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.0000000033025 0.00005747
## subject (Intercept) 0.0000000081275 0.00009015
## Residual 0.0000000006357 0.00002521
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.00124259 0.00002886 43.05
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.274
## subject | 0.674
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98089, p-value = 0.2759
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93303, p-value = 0.4423
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.91867, p-value = 0.006994
##
##
## $qq.resid
## [1] 18 17
##
## $qq.ranef.subject
## [1] 6 9
##
## $`qq.ranef.site:subject`
## [1] 17 12
##
## $res.v.fit
perf_mixed_model(models$ATROPOS_WM_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_WM_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1569.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.73497 -0.46998 -0.01091 0.49515 1.76890
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000002355 0.000004853
## subject (Intercept) 0.00000000040997 0.000020248
## Residual 0.00000000000382 0.000001954
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.00047857 0.00000616 77.69
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.054
## subject | 0.937
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.9918, p-value = 0.9183
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.91991, p-value = 0.3179
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.93204, p-value = 0.01878
##
##
## $qq.resid
## 38 37
## 37 36
##
## $qq.ranef.subject
## [1] 1 11
##
## $`qq.ranef.site:subject`
## [1] 21 30
##
## $res.v.fit
perf_mixed_model(models$ATROPOS_GM_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: ATROPOS_GM_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1502.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.64340 -0.40411 0.02559 0.35378 2.07830
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000006018 0.000007757
## subject (Intercept) 0.00000000024306 0.000015590
## Residual 0.00000000001381 0.000003716
## Number of obs: 74, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000693166 0.000004886 141.9
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.190
## subject | 0.767
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.92446, p-value = 0.0002792
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.9163, p-value = 0.289
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.98143, p-value = 0.7422
##
##
## $qq.resid
## 14 51
## 14 49
##
## $qq.ranef.subject
## [1] 8 1
##
## $`qq.ranef.site:subject`
## [1] 32 12
##
## $res.v.fit
perf_mixed_model(models$FAST_WM_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_WM_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1675.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.78643 -0.43514 0.02686 0.39648 1.86959
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000000000038492 0.000006204
## subject (Intercept) 0.000000000492182 0.000022185
## Residual 0.000000000005604 0.000002367
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000468050 0.000006769 69.14
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.072
## subject | 0.918
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.99336, p-value = 0.9574
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.91766, p-value = 0.2996
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.95814, p-value = 0.1447
##
##
## $qq.resid
## [1] 38 37
##
## $qq.ranef.subject
## [1] 1 11
##
## $`qq.ranef.site:subject`
## [1] 21 18
##
## $res.v.fit
perf_mixed_model(models$FAST_GM_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FAST_GM_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1642.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6925 -0.3926 0.0002 0.3444 2.4629
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000005085 0.000007131
## subject (Intercept) 0.00000000021045 0.000014507
## Residual 0.00000000001245 0.000003528
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000671595 0.000004541 147.9
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.186
## subject | 0.769
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.90592, p-value = 0.00002186
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.95301, p-value = 0.6826
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.97631, p-value = 0.555
##
##
## $qq.resid
## [1] 51 14
##
## $qq.ranef.subject
## [1] 1 8
##
## $`qq.ranef.site:subject`
## [1] 32 18
##
## $res.v.fit
perf_mixed_model(models$JLF_WM_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_WM_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1682.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.12628 -0.51321 -0.01984 0.35191 2.97520
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000000000019247 0.000004387
## subject (Intercept) 0.000000000469674 0.000021672
## Residual 0.000000000007354 0.000002712
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.00048768 0.00000658 74.12
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.039
## subject | 0.946
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.9634, p-value = 0.02177
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.9182, p-value = 0.3039
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.91105, p-value = 0.004072
##
##
## $qq.resid
## [1] 35 36
##
## $qq.ranef.subject
## [1] 1 11
##
## $`qq.ranef.site:subject`
## [1] 22 28
##
## $res.v.fit
perf_mixed_model(models$JLF_GM_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_GM_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1636.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.54869 -0.36475 0.00208 0.50417 1.80211
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.000000000093489 0.000009669
## subject (Intercept) 0.000000000517566 0.000022750
## Residual 0.000000000007899 0.000002811
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000744510 0.000007042 105.7
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.151
## subject | 0.836
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98258, p-value = 0.3476
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.93028, p-value = 0.4137
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.98775, p-value = 0.9365
##
##
## $qq.resid
## [1] 58 57
##
## $qq.ranef.subject
## [1] 8 1
##
## $`qq.ranef.site:subject`
## [1] 6 18
##
## $res.v.fit
perf_mixed_model(models$FIRST_THALAMUS_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: FIRST_THALAMUS_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1600.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.64679 -0.48633 -0.05519 0.53507 2.06990
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.00000000005445 0.000007379
## subject (Intercept) 0.00000000022682 0.000015060
## Residual 0.00000000002980 0.000005459
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000554717 0.000004736 117.1
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.175
## subject | 0.729
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.98806, p-value = 0.6687
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.98161, p-value = 0.9745
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.93314, p-value = 0.02042
##
##
## $qq.resid
## [1] 35 21
##
## $qq.ranef.subject
## [1] 7 2
##
## $`qq.ranef.site:subject`
## [1] 22 29
##
## $res.v.fit
perf_mixed_model(models$JLF_THALAMUS_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: JLF_THALAMUS_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1564.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1758 -0.3921 0.0092 0.4768 2.4299
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.0000000001390 0.000011792
## subject (Intercept) 0.0000000002651 0.000016283
## Residual 0.0000000000382 0.000006181
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.000539595 0.000005311 101.6
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.314
## subject | 0.599
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.94219, p-value = 0.001264
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.81838, p-value = 0.01642
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.97448, p-value = 0.4927
##
##
## $qq.resid
## [1] 51 52
##
## $qq.ranef.subject
## [1] 1 6
##
## $`qq.ranef.site:subject`
## [1] 22 10
##
## $res.v.fit
perf_mixed_model(models$MIMOSA_LESION_RD)
## $summary
## Linear mixed model fit by REML ['lmerMod']
## Formula: MIMOSA_LESION_RD ~ (1 | subject) + (1 | site:subject)
## Data: df
##
## REML criterion at convergence: -1352.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.36040 -0.45821 -0.06635 0.43847 1.65410
##
## Random effects:
## Groups Name Variance Std.Dev.
## site:subject (Intercept) 0.0000000028530 0.00005341
## subject (Intercept) 0.0000000031759 0.00005635
## Residual 0.0000000004712 0.00002171
## Number of obs: 80, groups: site:subject, 40; subject, 11
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.0007483 0.0000192 38.98
##
## $perf
## # ICC by Group
##
## Group | ICC
## --------------------
## site:subject | 0.439
## subject | 0.489
##
## $shapiro.resid
##
## Shapiro-Wilk normality test
##
## data: residuals(model)
## W = 0.97547, p-value = 0.1267
##
##
## $shapiro.subject
##
## Shapiro-Wilk normality test
##
## data: coef(model)$subject[, 1]
## W = 0.9102, p-value = 0.2452
##
##
## $`shapiro.site:subject`
##
## Shapiro-Wilk normality test
##
## data: coef(model)$`site:subject`[, 1]
## W = 0.99057, p-value = 0.9811
##
##
## $qq.resid
## [1] 24 19
##
## $qq.ranef.subject
## [1] 11 4
##
## $`qq.ranef.site:subject`
## [1] 14 33
##
## $res.v.fit
library(rstatix)
find_icc <- function(data){
irr::icc(data,
model = "twoway",
type = "agreement",
unit = "single")
}
run_anova <- function (var){
df %>%
filter(scan == '01') %>%
anova_test(dv = var, wid = subject, within = site) %>%
get_anova_table()
}
anovas <- outcomes %>%
lapply(run_anova) %>%
setNames(outcomes)
data.frame(measure = names(anovas),
F_stat = anovas %>%
purrr::map(~ .x$F) %>%
purrr::reduce(c),
p = anovas %>%
purrr::map(~ .x$p) %>%
purrr::reduce(c))
pair_df_list <- df %>%
pivot_longer(where(is.numeric)) %>%
split(list(.$site, .$name))
# compute names for each pair_df
df_names <- pair_df_list %>%
purrr::map(~ .x %>%
transmute(df_name = paste(site, name, sep="-")) %>%
pull() %>%
unique()
) %>%
purrr::reduce(c)
# finish transformation
pair_df_list <- pair_df_list %>%
purrr::map(~ .x %>%
pivot_wider(names_from = c('name', 'scan')) %>%
select(where(is.numeric))
) %>%
setNames(df_names)
icc_df <- pair_df_list %>%
purrr::map(find_icc) %>%
purrr::imap(~ data.frame(outcome = .y,
icc = .x$value,
lwr = .x$lbound,
upr = .x$ubound,
p_val = .x$p.value)) %>%
dplyr::bind_rows(.id = 'outcome') %>%
separate(outcome, into = c('site', 'outcome'), sep='-')
icc_df